******************************************************************************************************************************************
*** REPLICATION FILE FOR "ON THE ECONOMIC CONSEQUENCES OF MASS SHOOTINGS" BY ABEL BRODEUR AND HASIN YOUSAF
*** PREPARED FOR PUBLICATION IN THE REVIEW OF ECONOMICS AND STATISTICS
******************************************************************************************************************************************

clear all
set more off
set mat 3000
global path = "Put your path here"

******************************************************************************************************************************************
*** TABLE 1 
******************************************************************************************************************************************

*** - FIRST TWO COLUMNS 
use "$path\finaldata.dta", clear
sort fips year

foreach var of varlist ln_ca1_pop_1 ln_emp_pop ln_real_qp1_pop ln_est_pop ln_real_qp1_job ln_job_est lnhpi pop25g white black hispanic lhs coll {
reg L1.`var' mse ca1_pop_1 if aroundms==1 & postattack==0
}

foreach var of varlist ln_ca1_pop_1 ln_emp_pop ln_real_qp1_pop ln_est_pop ln_real_qp1_job ln_job_est lnhpi pop25g  white black hispanic lhs coll {
g ch`var'=D.`var'
}

foreach var of varlist ln_ca1_pop_1 ln_emp_pop ln_real_qp1_pop ln_est_pop ln_real_qp1_job ln_job_est lnhpi pop25g white black hispanic lhs  coll {
reg L1.ch`var' mse ca1_pop_1 if aroundms==1 & postattack==0
}

*** - LAST COLUMN
use "$path\finaldata.dta", clear

sort fips year

summ ln_ca1_pop_1 ln_emp_pop ln_real_qp1_pop ln_est_pop ln_real_qp1_job ln_job_est lnhpi pop25g white black hispanic lhs coll if aroundms==1 & postattack==0

******************************************************************************************************************************************
*** TABLE 2 
******************************************************************************************************************************************

use "$path\finaldata.dta", clear

reghdfe ln_emp_pop successful if aroundms==1, absorb(fips year) cluster(fips)
reghdfe ln_emp_pop successful if aroundms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_emp_pop successful if aroundms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_real_qp1_pop successful if aroundms==1, absorb(fips year) cluster(fips)
reghdfe ln_real_qp1_pop successful if aroundms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_real_qp1_pop successful if aroundms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_est_pop successful if aroundms==1, absorb(fips year) cluster(fips)
reghdfe ln_est_pop successful if aroundms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_est_pop successful if aroundms==1, absorb(fips i.div_9_all##i.year) cluster(fips)


******************************************************************************************************************************************
*** TABLE 3 
******************************************************************************************************************************************

*** - PANEL A
use "$path\finaldata.dta", clear

reghdfe ln_emp_pop successful postattack if aroundsuccfailms==1, absorb(fips i.year) cluster(fips)
reghdfe ln_emp_pop successful postattack if aroundsuccfailms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_emp_pop successful postattack if aroundsuccfailms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_real_qp1_pop successful postattack if aroundsuccfailms==1, absorb(fips i.year) cluster(fips)
reghdfe ln_real_qp1_pop successful postattack if aroundsuccfailms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_real_qp1_pop successful postattack if aroundsuccfailms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_est_pop successful postattack if aroundsuccfailms==1, absorb(fips i.year) cluster(fips)
reghdfe ln_est_pop successful postattack if aroundsuccfailms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_est_pop successful postattack if aroundsuccfailms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_emp_pop successful postattack if aroundsuccfailms==1, absorb(fips i.state_all##i.year) cluster(fips)
reghdfe ln_real_qp1_pop successful postattack if aroundsuccfailms==1, absorb(fips i.state_all##i.year) cluster(fips)
reghdfe ln_est_pop successful postattack if aroundsuccfailms==1, absorb(fips i.state_all##i.year) cluster(fips)


*** - PANEL B
use "$path\finaldata.dta", clear

reghdfe ln_emp_pop successful if aroundmsneighbor==1, absorb(fips year) cluster(fips)
reghdfe ln_emp_pop successful if aroundmsneighbor==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_emp_pop successful if aroundmsneighbor==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_real_qp1_pop successful if aroundmsneighbor==1, absorb(fips year) cluster(fips)
reghdfe ln_real_qp1_pop successful if aroundmsneighbor==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_real_qp1_pop successful if aroundmsneighbor==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_est_pop successful if aroundmsneighbor==1, absorb(fips year) cluster(fips)
reghdfe ln_est_pop successful if aroundmsneighbor==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_est_pop successful if aroundmsneighbor==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_emp_pop successful if aroundmsneighbor==1, absorb(fips i.state_all##i.year) cluster(fips)
reghdfe ln_real_qp1_pop successful if aroundmsneighbor==1, absorb(fips i.state_all##i.year) cluster(fips)
reghdfe ln_est_pop successful if aroundmsneighbor==1, absorb(fips i.state_all##i.year) cluster(fips)

*** - PANEL C
use "$path\FinalData\finaldata.dta", clear

sort fips year
foreach var of varlist ln_emp_pop ln_real_qp1_pop ln_est_pop ln_ca1_pop_1 lnVIOLENT lnPROPERTY male hs racediv poverty gini lnhom_gun lnsui_gun {
g lag`var'=L1.`var'
}

global controls = "lagln_emp_pop lagln_real_qp1_pop lagln_est_pop laglnVIOLENT laglnPROPERTY lagmale laghs lagracediv lagpoverty laggini laglnhom_gun laglnsui_gun capital_state coastal_county large_hub medium_hub"

***
probit ms lagln_ca1_pop_1 $controls i.year, cluster(fips)
predict p_score, pr
sort p_score year lagln_emp_pop lagln_real_qp1_pop lagln_est_pop fips

g matchms = 0
replace matchms = 1 if ms[_n+1]==1 & mse!=1 
replace matchms = 1 if ms[_n+2]==1 & mse!=1 
replace matchms = 1 if ms[_n+3]==1 & mse!=1 

sort fips year
by fips: egen matchmse = max(matchms)

g aroundmatchms = 0
forvalues j=-6(1)6{
by fips: replace aroundmatchms=1 if ms[_n+`j']==1
by fips: replace aroundmatchms=1 if matchms[_n+`j']==1
}

reghdfe ln_emp_pop successful if aroundmatchms==1, absorb(fips year) cluster(fips)
reghdfe ln_emp_pop successful if aroundmatchms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_emp_pop successful if aroundmatchms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_real_qp1_pop successful if aroundmatchms==1, absorb(fips year) cluster(fips)
reghdfe ln_real_qp1_pop successful if aroundmatchms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_real_qp1_pop successful if aroundmatchms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_est_pop successful if aroundmatchms==1, absorb(fips year) cluster(fips)
reghdfe ln_est_pop successful if aroundmatchms==1, absorb(fips i.region_4_all##i.year) cluster(fips)
reghdfe ln_est_pop successful if aroundmatchms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

reghdfe ln_emp_pop successful if aroundmatchms==1, absorb(fips i.state_all##i.year) cluster(fips)
reghdfe ln_real_qp1_pop successful if aroundmatchms==1, absorb(fips i.state_all##i.year) cluster(fips)
reghdfe ln_est_pop successful if aroundmatchms==1, absorb(fips i.state_all##i.year) cluster(fips)

************************************************************************************************************************
*** TABLE 4
************************************************************************************************************************

*** PANEL A
*Data not publicly available
*use "$path\sentiments.dta", clear
*global indcontrols "age age_2 male elem_school atthigh_school att_coll college married widowed divorced"

*reghdfe pers_fin_worse successful $indcontrols [pw=wt] if sample24==1, absorb(fips i.div_9_all##i.t) cluster(fips)
*reghdfe business_worse successful $indcontrols [pw=wt] if sample24==1, absorb(fips i.div_9_all##i.t) cluster(fips)
*reghdfe bad_time_hhitems successful $indcontrols [pw=wt] if sample24==1, absorb(fips i.div_9_all##i.t) cluster(fips)
*reghdfe pers_fin_worse_ft successful $indcontrols [pw=wt] if sample24==1, absorb(fips i.div_9_all##i.t) cluster(fips)

*** PANEL B
use "$path\BRFSS.dta", clear

reghdfe poorhlth successful elem_school atthigh_school high_school att_coll married divorced widowed separated couple male age age_2 [pw=_finalwt] ///
if sample==1 & retired!=1 & student!=1 & disabled!=1 & homemaker!=1  & age<=65 & menthlth!=., cluster(fips) absorb(fips year)

reghdfe physhlth successful elem_school atthigh_school high_school att_coll married divorced widowed separated couple male age age_2 [pw=_finalwt] ///
if sample==1 & retired!=1 & student!=1 & disabled!=1 & homemaker!=1  & age<=65 & menthlth!=., cluster(fips) absorb(fips year)

reghdfe menthlth successful elem_school atthigh_school high_school att_coll married divorced widowed separated couple male age age_2 [pw=_finalwt] ///
if sample==1 & retired!=1 & student!=1 & disabled!=1 & homemaker!=1 & age<=65, cluster(fips) absorb(fips year)

use "$path\finaldata.dta", clear
reghdfe lnhpi successful if aroundms==1, absorb(fips i.div_9_all##i.year) cluster(fips)

************************************************************************************************************************
*** TABLE 5
************************************************************************************************************************
use "$path\media_data2.dta", clear

*** PANEL A
***OLS
xi: reg chln_emp_pop media Total_city_minus_ms i.state_numeric i.division i.year succdeaths, cluster(fips) // .081

***IV
xi: ivreg2 chln_emp_pop (media = nd2) Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) first 
weakivtest

xi: reg chln_real_qp1_pop media Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) // .081

xi: ivreg2 chln_real_qp1_pop (media = nd2) Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) first 
weakivtest

xi: reg chln_est_pop media Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) // .081

xi: ivreg2 chln_est_pop (media = nd2) Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) first 
weakivtest

*** PANEL B
***OLS
xi: reg chln_emp_pop media_min Total_city_minus_ms i.state_numeric i.division i.year succdeaths, cluster(fips) // .081

***IV
xi: ivreg2 chln_emp_pop (media_min = nd2) Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) first 
weakivtest

xi: reg chln_real_qp1_pop media_min Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) // .081

xi: ivreg2 chln_real_qp1_pop (media_min = nd2) Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) first 
weakivtest

xi: reg chln_est_pop media_min Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) // .081

xi: ivreg2 chln_est_pop (media_min = nd2) Total_city_minus_ms i.state_numeric i.divisionyear succdeaths, cluster(fips) first 
weakivtest

******************************************************************************************************************************************
*** FIGURES
******************************************************************************************************************************************

******************************************************************************************************************************************
*** FIGURE 1 - MAP
******************************************************************************************************************************************
clear
*Change path to a folder with all Shapefiles or use
*shp2dta using "$path\gz_2010_us_050_00_500k", data("$path\gz_2010_us_050_00_500k.dta")  coor("$path\gz_2010_us_050_00_500k_coor.dta") replace

shp2dta using "$path\Shapefile\gz_2010_us_050_00_500k", data("$path\Shapefile\gz_2010_us_050_00_500k.dta")  coor("$path\Shapefile\gz_2010_us_050_00_500k_coor.dta") replace

use "$path\Shapefile\gz_2010_us_050_00_500k.dta", clear
g fips=1000*real(STATE)+real(COUNTY)
compress
save "$path\Shapefile\gz_2010_us_050_00_500k.dta", replace

use "$path\finaldata.dta" , clear

sort fips year
foreach var of varlist ln_emp_pop ln_real_qp1_pop ln_est_pop ln_ca1_pop_1 lnVIOLENT lnPROPERTY male hs racediv poverty gini lnhom_gun lnsui_gun {
g lag`var'=L1.`var'
}

global controls = "lagln_emp_pop lagln_real_qp1_pop lagln_est_pop laglnVIOLENT laglnPROPERTY lagmale laghs lagracediv lagpoverty laggini laglnhom_gun laglnsui_gun capital_state coastal_county large_hub medium_hub"

***
probit ms lagln_ca1_pop_1 $controls i.year, cluster(fips)
predict p_score, pr
sort p_score year lagln_emp_pop lagln_real_qp1_pop lagln_est_pop fips

g matchms = 0
replace matchms = 1 if ms[_n+1]==1 & mse!=1 
replace matchms = 1 if ms[_n+2]==1 & mse!=1 
replace matchms = 1 if ms[_n+3]==1 & mse!=1 

sort fips year
by fips: egen matchmse = max(matchms)

duplicates drop fips mse failmse matchmse, force
keep fips mse failmse matchmse
duplicates drop fips, force
merge 1:1 fips using "$path\Shapefile\gz_2010_us_050_00_500k.dta"
g category = 1*(mse==1 & failmse==0 & matchmse==0) + 2*(mse==0 & failmse==1 & matchmse==0) + 3*(mse==0 & failmse==0 & matchmse==1) + 4*(mse==1 & failmse==1 & matchmse==0) + 5*(mse==0 & failmse==1 & matchmse==1)
replace category=. if category==0
label define category 1 "MS Only" 2 "Failed MS" 3 "Matched Counties" 4 "MS & Failed MS" 5 "Failed MS & Matched Counties" , modify
label values category category 
drop if STATE=="02" | STATE=="15" | STATE=="72"

g category2 = 1*(mse==1 & failmse==0 & matchmse==0) + 2*(mse==0 & (failmse==1 | matchmse==1)) + 3*(mse==1 & (failmse==1 | matchmse==1))
label define category2 0 "None" 1 "MS Only" 2 "Failed or Matched MS" 3 "MS & Failed or Matched MS"
label values category2 category2 

spmap category using "$path\Shapefile\gz_2010_us_050_00_500k_coor.dta", id(_ID) clmethod(unique) legtitle("Mass Shooting") leglabel("category") fcolor(RdYlGn)  

******************************************************************************************************************************************
*** FIGURE 2 
******************************************************************************************************************************************
use "$path\finaldata.dta", clear
sort fips year

g beta_succ_emp = .
g se_succ_emp = .

reghdfe ln_emp_pop successful13-successful17 successful19-successful25 if aroundms==1, absorb(fips i.div_9_all##i.year) cluster(fips)
foreach var of varlist successful13-successful17 successful19-successful25 {
replace beta_succ_emp = _b[`var'] if `var'==1
replace se_succ_emp = _se[`var'] if `var'==1
}
replace beta_succ_emp = 0 if successful18==1
replace se_succ_emp = 0 if successful18==1
g lbar_succ_emp = beta_succ_emp - 1.64*se_succ_emp
g ubar_succ_emp = beta_succ_emp + 1.64*se_succ_emp

g beta_succ_ern = .
g se_succ_ern = .

reghdfe ln_real_qp1_pop successful13-successful17 successful19-successful25 if aroundms==1, absorb(fips i.div_9_all##i.year) cluster(fips)
foreach var of varlist successful13-successful17 successful19-successful25 {
replace beta_succ_ern = _b[`var'] if `var'==1
replace se_succ_ern = _se[`var'] if `var'==1
}
replace beta_succ_ern = 0 if successful18==1
replace se_succ_ern = 0 if successful18==1
g lbar_succ_ern = beta_succ_ern - 1.64*se_succ_ern
g ubar_succ_ern = beta_succ_ern + 1.64*se_succ_ern

g beta_succ_est = .
g se_succ_est = .

reghdfe ln_est_pop successful13-successful17 successful19-successful25 if aroundms==1, absorb(fips i.div_9_all##i.year) cluster(fips)
foreach var of varlist successful13-successful17 successful19-successful25 {
replace beta_succ_est = _b[`var'] if `var'==1
replace se_succ_est = _se[`var'] if `var'==1
}
replace beta_succ_est = 0 if successful18==1
replace se_succ_est = 0 if successful18==1
g lbar_succ_est = beta_succ_est - 1.64*se_succ_est
g ubar_succ_est = beta_succ_est + 1.64*se_succ_est

keep if yearsince>=-5 & yearsince<=5
duplicates drop yearsince, force
foreach var of varlist lbar* ubar* {
replace `var'=0 if `var'==.
}

twoway (connected beta_succ_emp yearsince) (line lbar_succ_emp yearsince, lpattern(dash) lcolor(black)) ///
(line ubar_succ_emp yearsince, lpattern(dash) lcolor(black) yline(0, lcolor(black))), ///
xtitle("Years since Mass Shooting") ytitle("Impact of Successful MS on Employment") ///
ylabel(-5(1)2) xlabel(-5(1)5) legend(off)

twoway (connected beta_succ_ern yearsince) (line lbar_succ_ern yearsince, lpattern(dash) lcolor(black)) ///
(line ubar_succ_ern yearsince, lpattern(dash) lcolor(black) yline(0, lcolor(black))), ///
xtitle("Years since Mass Shooting") ytitle("Impact of Successful MS on Earnings") ///
ylabel(-5(1)2) xlabel(-5(1)5) legend(off)

twoway (connected beta_succ_est yearsince) (line lbar_succ_est yearsince, lpattern(dash) lcolor(black)) ///
(line ubar_succ_est yearsince, lpattern(dash) lcolor(black) yline(0, lcolor(black))), ///
xtitle("Years since Mass Shooting") ytitle("Impact of Successful MS on Establishments") ///
ylabel(-5(1)2) xlabel(-5(1)4) legend(off)
